{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "import os,sys,inspect\n",
    "import os\n",
    "import joblib\n",
    "import tensorflow as tf\n",
    "import numpy as np\n",
    "import h5py\n",
    "import scipy.sparse.linalg as la\n",
    "import scipy.sparse as sp\n",
    "import scipy\n",
    "import time\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "from matplotlib.backends.backend_pdf import PdfPages\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "path_dataset = '../../datasets/douban/training_test_dataset.mat' #with only social network"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "# auxiliary functions:\n",
    "\n",
    "# import matlab files in python\n",
    "def load_matlab_file(path_file, name_field):\n",
    "    \"\"\"\n",
    "    load '.mat' files\n",
    "    inputs:\n",
    "        path_file, string containing the file path\n",
    "        name_field, string containig the field name (default='shape')\n",
    "    warning:\n",
    "        '.mat' files should be saved in the '-v7.3' format\n",
    "    \"\"\"\n",
    "    db = h5py.File(path_file, 'r')\n",
    "    ds = db[name_field]\n",
    "    try:\n",
    "        if 'ir' in ds.keys():\n",
    "            data = np.asarray(ds['data'])\n",
    "            ir   = np.asarray(ds['ir'])\n",
    "            jc   = np.asarray(ds['jc'])\n",
    "            out  = sp.csc_matrix((data, ir, jc)).astype(np.float32)\n",
    "    except AttributeError:\n",
    "        # Transpose in case is a dense matrix because of the row- vs column- major ordering between python and matlab\n",
    "        out = np.asarray(ds).astype(np.float32).T\n",
    "\n",
    "    db.close()\n",
    "\n",
    "    return out"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "#loading of the required matrices\n",
    "M = load_matlab_file(path_dataset, 'M')\n",
    "Otraining = load_matlab_file(path_dataset, 'Otraining')\n",
    "Otest = load_matlab_file(path_dataset, 'Otest')\n",
    "Wrow = load_matlab_file(path_dataset, 'W_users') #dense"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "np.random.seed(0)\n",
    "pos_tr_samples = np.where(Otraining)\n",
    "\n",
    "num_tr_samples = len(pos_tr_samples[0])\n",
    "list_idx = range(num_tr_samples)\n",
    "np.random.shuffle(list_idx)\n",
    "idx_data = list_idx[:num_tr_samples//2]\n",
    "idx_train = list_idx[num_tr_samples//2:]\n",
    "\n",
    "pos_data_samples = (pos_tr_samples[0][idx_data], pos_tr_samples[1][idx_data])\n",
    "pos_tr_samples = (pos_tr_samples[0][idx_train], pos_tr_samples[1][idx_train])\n",
    "\n",
    "Odata = np.zeros(M.shape)\n",
    "Otraining = np.zeros(M.shape)\n",
    "\n",
    "for k in range(len(pos_data_samples[0])):\n",
    "    Odata[pos_data_samples[0][k], pos_data_samples[1][k]] = 1\n",
    "    \n",
    "for k in range(len(pos_tr_samples[0])):\n",
    "    Otraining[pos_tr_samples[0][k], pos_tr_samples[1][k]] = 1\n",
    "    \n",
    "print 'Num data samples: %d' % (np.sum(Odata),)\n",
    "print 'Num train samples: %d' % (np.sum(Otraining),)\n",
    "print 'Num train+data samples: %d' % (np.sum(Odata+Otraining),)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "#computation of the normalized laplacians\n",
    "Lrow = sp.csgraph.laplacian(Wrow, normed=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "#apply SVD initially for detecting the main components of our initialization\n",
    "U, s, V = np.linalg.svd(Odata*M, full_matrices=0)\n",
    "\n",
    "print U.shape\n",
    "print s.shape\n",
    "print V.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "rank_W_H = 10\n",
    "partial_s = s[:rank_W_H]\n",
    "partial_S_sqrt = np.diag(np.sqrt(partial_s))\n",
    "initial_W = np.dot(U[:, :rank_W_H], partial_S_sqrt)\n",
    "initial_H = np.dot(partial_S_sqrt, V[:rank_W_H, :]).T\n",
    "\n",
    "print initial_W.shape\n",
    "print initial_H.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "print 'Original training matrix'\n",
    "plt.figure()\n",
    "plt.imshow(Odata*M)\n",
    "plt.colorbar()\n",
    "\n",
    "print 'Reconstructed training matrix'\n",
    "plt.figure()\n",
    "plt.imshow(np.dot(initial_W, initial_H.T))\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "class Train_test_matrix_completion:\n",
    "    \n",
    "    \"\"\"\n",
    "    The neural network model.\n",
    "    \"\"\"\n",
    "    \n",
    "    def frobenius_norm(self, tensor):\n",
    "        square_tensor = tf.square(tensor)\n",
    "        tensor_sum = tf.reduce_sum(square_tensor)\n",
    "        frobenius_norm = tf.sqrt(tensor_sum)\n",
    "        return frobenius_norm\n",
    "    \n",
    "    def mono_conv(self, list_lap, ord_conv, A, W, b):\n",
    "        \n",
    "        feat = []\n",
    "        #collect features\n",
    "        for k in range(ord_conv):\n",
    "            c_lap = list_lap[k] \n",
    "                                                     \n",
    "            #dense implementation\n",
    "            c_feat = tf.matmul(c_lap, A, a_is_sparse=False)\n",
    "            feat.append(c_feat)\n",
    "            \n",
    "        all_feat = tf.concat(feat,1)\n",
    "        conv_feat = tf.matmul(all_feat, W) + b\n",
    "        conv_feat = tf.nn.relu(conv_feat)\n",
    "        \n",
    "        return conv_feat\n",
    "                \n",
    "    def compute_cheb_polynomials(self, L, ord_cheb, list_cheb):\n",
    "        for k in range(ord_cheb):\n",
    "            if (k==0):\n",
    "                list_cheb.append(tf.cast(tf.diag(tf.ones([tf.shape(L)[0],])), 'float32'))\n",
    "            elif (k==1):\n",
    "                list_cheb.append(tf.cast(L, 'float32'))\n",
    "            else:\n",
    "                list_cheb.append(2*tf.matmul(L, list_cheb[k-1])  - list_cheb[k-2])\n",
    "                 \n",
    "    def __init__(self, M, Lr, Odata, Otraining, Otest, initial_W, initial_H,\n",
    "                 order_chebyshev_row = 5,\n",
    "                 num_iterations = 10, gamma=1.0, gamma_H=1.0, gamma_W=1.0, learning_rate=1e-4, idx_gpu = '/gpu:1'):\n",
    "        \n",
    "        #order of the spectral filters\n",
    "        self.ord_row = order_chebyshev_row\n",
    "        self.num_iterations = num_iterations\n",
    "        self.n_conv_feat = 32\n",
    "        \n",
    "        with tf.Graph().as_default() as g:\n",
    "                tf.logging.set_verbosity(tf.logging.ERROR)\n",
    "                self.graph = g\n",
    "                tf.set_random_seed(0)\n",
    "                with tf.device(idx_gpu):\n",
    "                    \n",
    "                        #loading of the laplacians\n",
    "                        self.Lr = tf.cast(Lr, 'float32')\n",
    "                        \n",
    "                        self.norm_Lr = self.Lr - tf.diag(tf.ones([Lr.shape[0], ]))\n",
    "                        \n",
    "                        #compute all chebyshev polynomials a priori\n",
    "                        self.list_row_cheb_pol = list()\n",
    "                        self.compute_cheb_polynomials(self.norm_Lr, self.ord_row, self.list_row_cheb_pol)\n",
    "                        \n",
    "                        #definition of constant matrices\n",
    "                        self.M = tf.constant(M, dtype=tf.float32)\n",
    "                        self.Odata = tf.constant(Odata, dtype=tf.float32)\n",
    "                        self.Otraining = tf.constant(Otraining, dtype=tf.float32) #training mask\n",
    "                        self.Otest = tf.constant(Otest, dtype=tf.float32) #test mask\n",
    "                         \n",
    "                        ##################################definition of the NN variables#####################################\n",
    "                        \n",
    "                        #definition of the weights for extracting the global features\n",
    "                        self.W_conv_W = tf.get_variable(\"W_conv_W\", shape=[self.ord_row*initial_W.shape[1], self.n_conv_feat], initializer=tf.contrib.layers.xavier_initializer())\n",
    "                        self.b_conv_W = tf.Variable(tf.zeros([self.n_conv_feat,]))\n",
    "                        \n",
    "                        #recurrent N parameters\n",
    "                        self.W_f_u = tf.get_variable(\"W_f_u\", shape=[self.n_conv_feat, self.n_conv_feat], initializer=tf.contrib.layers.xavier_initializer())\n",
    "                        self.W_i_u = tf.get_variable(\"W_i_u\", shape=[self.n_conv_feat, self.n_conv_feat], initializer=tf.contrib.layers.xavier_initializer())\n",
    "                        self.W_o_u = tf.get_variable(\"W_o_u\", shape=[self.n_conv_feat, self.n_conv_feat], initializer=tf.contrib.layers.xavier_initializer())\n",
    "                        self.W_c_u = tf.get_variable(\"W_c_u\", shape=[self.n_conv_feat, self.n_conv_feat], initializer=tf.contrib.layers.xavier_initializer())\n",
    "                        self.U_f_u = tf.get_variable(\"U_f_u\", shape=[self.n_conv_feat, self.n_conv_feat], initializer=tf.contrib.layers.xavier_initializer())\n",
    "                        self.U_i_u = tf.get_variable(\"U_i_u\", shape=[self.n_conv_feat, self.n_conv_feat], initializer=tf.contrib.layers.xavier_initializer())\n",
    "                        self.U_o_u = tf.get_variable(\"U_o_u\", shape=[self.n_conv_feat, self.n_conv_feat], initializer=tf.contrib.layers.xavier_initializer())\n",
    "                        self.U_c_u = tf.get_variable(\"U_c_u\", shape=[self.n_conv_feat, self.n_conv_feat], initializer=tf.contrib.layers.xavier_initializer())\n",
    "                        self.b_f_u = tf.Variable(tf.zeros([self.n_conv_feat,]))\n",
    "                        self.b_i_u = tf.Variable(tf.zeros([self.n_conv_feat,]))\n",
    "                        self.b_o_u = tf.Variable(tf.zeros([self.n_conv_feat,]))\n",
    "                        self.b_c_u = tf.Variable(tf.zeros([self.n_conv_feat,]))\n",
    "                        \n",
    "                        #output parameters\n",
    "                        self.W_out_W = tf.get_variable(\"W_out_W\", shape=[self.n_conv_feat, initial_W.shape[1]], initializer=tf.contrib.layers.xavier_initializer()) \n",
    "                        self.b_out_W = tf.Variable(tf.zeros([initial_W.shape[1],]))\n",
    "                        \n",
    "                        #########definition of the NN\n",
    "                        #definition of W and H\n",
    "                        self.W = tf.constant(initial_W.astype('float32'))\n",
    "                        self.H = tf.Variable(initial_H.astype('float32'))\n",
    "                        \n",
    "                        self.X = tf.matmul(self.W, self.H, transpose_b=True) #we may initialize it at random here\n",
    "                        self.list_X = list()\n",
    "                        self.list_X.append(tf.identity(self.X))\n",
    "                        \n",
    "                        #RNN\n",
    "                        self.h_u = tf.zeros([M.shape[0], self.n_conv_feat])\n",
    "                        self.c_u = tf.zeros([M.shape[0], self.n_conv_feat])\n",
    "                        \n",
    "                        \n",
    "                        for k in range(self.num_iterations):\n",
    "                            #extraction of global features vectors\n",
    "                            self.final_feat_users = self.mono_conv(self.list_row_cheb_pol, self.ord_row, self.W, self.W_conv_W, self.b_conv_W)\n",
    "                            \n",
    "                            #here we have to split the features between users and movies LSTMs\n",
    "                            \n",
    "                            #users RNN\n",
    "                            self.f_u = tf.sigmoid(tf.matmul(self.final_feat_users, self.W_f_u) + tf.matmul(self.h_u, self.U_f_u) + self.b_f_u)\n",
    "                            self.i_u = tf.sigmoid(tf.matmul(self.final_feat_users, self.W_i_u) + tf.matmul(self.h_u, self.U_i_u) + self.b_i_u)\n",
    "                            self.o_u = tf.sigmoid(tf.matmul(self.final_feat_users, self.W_o_u) + tf.matmul(self.h_u, self.U_o_u) + self.b_o_u)\n",
    "                            \n",
    "                            self.update_c_u = tf.sigmoid(tf.matmul(self.final_feat_users, self.W_c_u) + tf.matmul(self.h_u, self.U_c_u) + self.b_c_u)\n",
    "                            self.c_u = tf.multiply(self.f_u, self.c_u) + tf.multiply(self.i_u, self.update_c_u)\n",
    "                            self.h_u = tf.multiply(self.o_u, tf.sigmoid(self.c_u))\n",
    "                            \n",
    "                            #compute update of matrix X\n",
    "                            self.delta_W = tf.tanh(tf.matmul(self.c_u, self.W_out_W) + self.b_out_W) #N x rank_W_H\n",
    "                            \n",
    "                            self.W += self.delta_W\n",
    "                        \n",
    "                            self.X = tf.matmul(self.W, self.H, transpose_b=True)\n",
    "                            self.list_X.append(tf.identity(tf.reshape(self.X, [tf.shape(self.M)[0], tf.shape(self.M)[1]])))\n",
    "                        self.X = tf.matmul(self.W, self.H, transpose_b=True)\n",
    "                        #########loss definition\n",
    "                        \n",
    "                        #computation of the accuracy term\n",
    "                        self.norm_X = 1+4*(self.X-tf.reduce_min(self.X))/(tf.reduce_max(self.X-tf.reduce_min(self.X)))\n",
    "                        frob_tensor = tf.multiply(self.Otraining + self.Odata, self.norm_X - M)\n",
    "                        self.loss_frob = tf.square(self.frobenius_norm(frob_tensor))/np.sum(Otraining+Odata)\n",
    "                        \n",
    "                        #computation of the regularization terms\n",
    "                        trace_row_tensor = tf.matmul(tf.matmul(self.X, self.Lr, transpose_a=True), self.X)\n",
    "                        self.loss_trace_row = tf.trace(trace_row_tensor)/tf.cast(tf.shape(self.X)[0]*tf.shape(self.X)[1],'float32')\n",
    "                        \n",
    "                        self.frob_norm_H = tf.square(self.frobenius_norm(self.H))/tf.cast(tf.shape(self.H)[0]*tf.shape(self.H)[1], 'float32')\n",
    "                        self.frob_norm_W = tf.square(self.frobenius_norm(self.W))/tf.cast(tf.shape(self.W)[0]*tf.shape(self.W)[1], 'float32')\n",
    "                        \n",
    "                        \n",
    "                        #training loss definition\n",
    "                        self.loss = self.loss_frob + (gamma/2)*self.loss_trace_row + (gamma_H/2)*self.frob_norm_H + (gamma_W/2)*self.frob_norm_W\n",
    "                        \n",
    "                        #test loss definition\n",
    "                        self.predictions = tf.multiply(self.Otest, self.norm_X - self.M)\n",
    "                        self.predictions_error = self.frobenius_norm(self.predictions)\n",
    "\n",
    "                        #definition of the solver\n",
    "                        self.optimizer = tf.train.AdamOptimizer(learning_rate=learning_rate).minimize(self.loss)\n",
    "                        \n",
    "                        self.var_grad = tf.gradients(self.loss, tf.trainable_variables())\n",
    "                        self.norm_grad = self.frobenius_norm(tf.concat([tf.reshape(g, [-1]) for g in self.var_grad],0))\n",
    "\n",
    "                        # Create a session for running Ops on the Graph.\n",
    "                        config = tf.ConfigProto(allow_soft_placement = True)\n",
    "                        config.gpu_options.allow_growth = True\n",
    "                        self.session = tf.Session(config=config)\n",
    "\n",
    "                        # Run the Op to initialize the variables.\n",
    "                        init = tf.initialize_all_variables()\n",
    "                        self.session.run(init)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "ord_col = 5\n",
    "ord_row = 5\n",
    "\n",
    "        \n",
    "learning_obj = Train_test_matrix_completion(M, Lrow, Odata, Otraining, Otest, \n",
    "                                                    initial_W, initial_H,\n",
    "                                                    order_chebyshev_row = ord_row, \n",
    "                                                    gamma=1e2, gamma_H=1e2, gamma_W=1e0,\n",
    "                                                    learning_rate=1e-3)\n",
    "\n",
    "num_iter_test = 10\n",
    "num_total_iter_training = 5000\n",
    "\n",
    "num_iter = 0\n",
    "\n",
    "list_training_loss = list()\n",
    "list_training_norm_grad = list()\n",
    "list_test_pred_error = list()\n",
    "list_predictions = list()\n",
    "list_X = list()\n",
    "\n",
    "list_training_times = list()\n",
    "list_test_times = list()\n",
    "list_grad_X = list()\n",
    "\n",
    "list_X_evolutions = list()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "num_iter = 0\n",
    "for k in range(num_iter, num_total_iter_training):\n",
    "\n",
    "    tic = time.time()\n",
    "    _, current_training_loss, norm_grad, X_grad = learning_obj.session.run([learning_obj.optimizer, learning_obj.loss, \n",
    "                                                                                        learning_obj.norm_grad, learning_obj.var_grad]) \n",
    "    training_time = time.time() - tic\n",
    "\n",
    "    list_training_loss.append(current_training_loss)\n",
    "    list_training_norm_grad.append(norm_grad)\n",
    "    list_training_times.append(training_time)\n",
    "\n",
    "    if (np.mod(num_iter, num_iter_test)==0):\n",
    "        msg = \"[TRN] iter = %03i, cost = %3.2e, |grad| = %.2e (%3.2es)\" \\\n",
    "                                    % (num_iter, list_training_loss[-1], list_training_norm_grad[-1], training_time)\n",
    "        print msg\n",
    "\n",
    "        #Test Code\n",
    "        tic = time.time()\n",
    "        pred_error, preds, X = learning_obj.session.run([learning_obj.predictions_error, learning_obj.predictions,\n",
    "                                                                             learning_obj.norm_X]) \n",
    "        test_time = time.time() - tic\n",
    "\n",
    "        list_test_pred_error.append(pred_error)\n",
    "        list_test_times.append(test_time)\n",
    "            \n",
    "        RMSE = np.sqrt(np.square(pred_error)/np.sum(Otest))\n",
    "        msg =  \"[TST] iter = %03i, cost = %3.2e, RMSE = %3.2e (%3.2es)\" % (num_iter, list_test_pred_error[-1], RMSE, test_time)\n",
    "        print msg\n",
    "            \n",
    "    num_iter += 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "fig, ax1 = plt.subplots(figsize=(20,10))\n",
    "\n",
    "ax2 = ax1.twinx()\n",
    "ax1.plot(np.arange(len(list_training_loss)), list_training_loss, 'g-')\n",
    "ax2.plot(np.arange(len(list_test_pred_error))*num_iter_test, list_test_pred_error, 'b-')\n",
    "\n",
    "ax1.set_xlabel('Iteration')\n",
    "ax1.set_ylabel('Training loss', color='g')\n",
    "ax2.set_ylabel('Test loss', color='b')\n",
    "\n",
    "best_iter = (np.where(np.asarray(list_training_loss)==np.min(list_training_loss))[0][0]//num_iter_test)*num_iter_test\n",
    "best_pred_error = list_test_pred_error[best_iter//num_iter_test]\n",
    "\n",
    "print 'Best predictions at iter: %d (error: %f)' % (best_iter, best_pred_error)\n",
    "RMSE = np.sqrt(np.square(best_pred_error)/np.sum(Otest))\n",
    "print 'RMSE: %f' % RMSE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "#last X generated\n",
    "plt.figure(figsize=(20,10))\n",
    "plt.imshow(X)\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
